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Abstract 

Q\ ' This paper presents comparison of several stochastic optimization algorithms 

developed by authors in their previous works for the solution of some problems 
^ ■ arising in Civil Engineering. The introduced optimization methods are: the inte- 

! ger augmented simulated anneahng (lASA), the real-coded augmented simulated 

^ I annealing (RASA) [TD] , the differential evolution (DE) in its original fashion de- 

veloped by R. Storn and K. Price [T3] and simplified real-coded differential genetic 



X 



algorithm (SADE) [6J. Each of these methods was developed for some specific 
optimization problem; namely the Chebychev trial polynomial problem, the so 
! called type function and two engineering problems - the reinforced concrete 

^ I beam layout and the periodic unit cell problem respectively. Detailed and ex- 

tyj I tensive numerical tests were performed to examine the stability and efficiency 

of proposed algorithms. The results of our experiments suggest that the per- 
formance and robustness of RASA, lASA and SADE methods are comparable, 
. while the DE algorithm performs slightly worse. This fact together with a small 

I number of internal parameters promotes the SADE method as the most robust 

^ I for practical use. 
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1 Introduction 

Nowadays, optimization has become one of the most discussed topics of engi- 
neering and applied research. Advantages coming from using optimization tools 
in engineering design are obvious. They allow to choose an optimal layout of 
a certain structure or a structural component from the huge space of possible 
solutions based on a more realistic physical model, while the traditional design- 
ing methods usually rely only on some simple empirical formulas or guidelines 
resulting in a feasible but not necessarily an (sub-)optimal solution. Using opti- 
mization as a method of design can raise engineering job to a higher level, both 
in terms of efficiency and reliability of obtained results. 



^Corresponding author, e-mail: lepsOcml . f sv . cvut . cz, fax: +420-2-2431-077 
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Typically, optimization methods arising in engineering design problems are 
computationally demanding because they require evaluation of a quite compli- 
cated objective function many times for different potential solutions. Moreover, 
the objective function is often multi-modal, non-smooth or even discontinuous, 
which means that traditional, gradient-based optimization algorithms fail and 
global optimization techniques, which generally need even a larger number of 
function calls, must be employed. Fortunately, the rapid development of compu- 
tational technologies and hardware components allows us to treat these problems 
within a reasonable time. 

As indicated previously, the optimization methods can be divided generally 
into two groups: the gradient methods, that operate on a single potential solution 
and look for some improvements in its neighborhood, and global optimization 
techniques - represented here by so called evolutionary methods - that maintain 
large sets (populations) of potential solutions and apply some recombination 
and selection operators on them. During the last decades, evolutionary methods 
have received a considerable attraction and have experienced a rapid develop- 
ment. Good compendium of these methods can be found for example in [12] and 
references therein. Main paradigms are: genetic algorithm (binary or real coded), 
augmented simulated annealing (binary or real coded), evolution strategy and dif- 
ferential evolution. Each of these methods has many possible improvements (see, 

e.g., mM)- 

Many researchers all over the world are united in an effort to develop an 
evolutionary optimization method that is able to solve reliably any problem 
submitted to it. In present time, there is no such method available. Each method 
is able to outperform others for certain type of optimization problem, but it 
extremely slows down or even fails for another one. Moreover, many authors do 
not introduce the reliable testing methodology for ranking their methods. For 
example they introduce results of a single run of a given method, which is rather 
questionable for the case of stochastic algorithms. Finally, the methods are often 
benchmarked on some test functions, that even if presented as complicated, are 
continuous and have few local extremes. 

This paper presents several optimization methods that were developed and 
tested for different types of optimization tasks. The integer augmented simulated 
annealing (lASA), derived from a binary version of the algorithm P], was devel- 
oped to optimize a reinforced concrete beam layout from the economic point of 
view. For solving the problem of a periodic unit cell layout [16j, the real-coded 
simulated annealing was applied. Differential evolution arose to solve famous 
Chebychev polynomial problem The last of the introduced methods is 

the so-called SADE technology. It is a simplified real-coded differential genetic 
algorithms that was developed as a specific recombination of a genetic algorithm 
and a differential evolution intended to solve problems on high-dimensional do- 
mains represented by the type Otest function |6j. All these methods may aspire to 
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be a universal optimization algorithm. So, we have conducted a detailed numer- 
ical tests of all these four optimization methods for aforementioned optimization 
problems to examine their behavior and performance. 

The paper is organized as follows. Section 2 provides brief description of 
each optimization task, while individual optimization algorithms are discussed in 
Section 3. Numerical results appear in Section 4. Summary on the performance 
of individual methods is presented in Section 5. For the sake of completeness, 
the parameter settings of algorithms is shown in the Appendix. 

2 Optimization tasks 

The optimization problems that are used as a set of test functions can be divided 
into two groups: the "test suite", containing "artificial functions" and the "engi- 
neering problems", which collect more (hopefully) practical optimization tasks. 
Specifically, these problems are : 

• Test suite: 

— Chebychev trial polynomial problem, 

— Type benchmark trial function. 

• Engineering problems 

— Reinforced concrete beam layout, 

— Periodic unit cell problem. 

The following section provides description of selected functions in more details. 
2.1 Chebychev problem 

Chebychev trial polynomial problem is one of the most famous optimization 
problems. Our goal is to find such coefficients of a polynomial constrained by 
the condition that the graph of the polynomial can be fitted into a specified area 
(see Fig. [ID . 

Thus, the optimized values are the parameters a, of a polynomial expression: 

n 

f{x)=Y,aix\ (1) 

and the value of objective function is determined as a sum of the areas, where 
the function graph exceeds a given boundary (hatched areas in Fig. [1]). 
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Figure 1: A graph of a Chebychev polynomial (n = 8). 



2.2 T)jpe function 



This trial optimization problem was proposed by the first two authors to examine 
the ability of the optimization method to find a single extreme of a function with 
a high number of parameters and growth of computational complexity with the 
problem dimension. For this reason, we used a function with a single extreme 
on the top of the high and narrow peak: 



where x is a vector of unknown variables, xq is the point of the global extreme 
(the top of the peak) and Uq and tq are parameters that infiuence the height or the 
width of the peak, respectively. Example of such a function on one dimensional 
domain is shown in Fig. [21 

Although this example function has only a single extreme, to find it even with 
a moderate precision is a non-trivial task for several reasons. First, in the very 
neighborhood of the extreme the function is so steep that even a futile change 
of the coordinates cause a large change of the function value; in such a case it is 
very difficult for the algorithm to determine what way leads to the top. Second, 
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Figure 2: An example of a type function. 



the peak is located on a very narrow part of a domain and this disproportion 
increases very quickly with the problem dimension. 

2.3 Reinforced concrete beam layout 

An effort to create an optimal design of a steel-reinforced concrete structure is 
as old as the material itself. In present times an emphasis is put on this problem 
due to widespread use of RC structures in Civil Engineering. Frame structures 
are major part in this field with beams playing an important role as one of 
the basic building block of this construction system. An objective is to choose 
the best design from all possible configurations that can create the requested 
structure - in our case a continuous beam (see Fig. [3l). 

The total cost of the structure is used as a comparison factor. An advantage 
of the financial rating is its natural meaning to non-experts and easiness of 
constraints implementation. In our particular case, the objective function reads 

f{X) = V,P, + W,P, + Y.pfi, (3) 
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Figure 3: A continuous beam subjected to an uniform loading. 

where is the volume of concrete and Ws is the weight of steel; Pc and Pg 
are the price of concrete per unit volume and steel per kilogram, respectively. 
From the mathematical point of view the penalty function pfi is a distance 
between a solution and the feasible space, or equivalently, a price spent on the 
fulfillment of the condition i. Suppose that a variable should not exceed a 
certain allowable limit $j,max- Then, the penalty pfi assumes the form 







if *i < ^i,max, 

otherwise, 



(4) 



where Wi is the weight of the i-th constraint. 

The constraints in this procedure deal with allowable strength and service- 
ability limits given by a chosen standard - in our case EUROCODE 2 (EC2) [2]. 
An interested reader can find implementation details for example in [0]. 

Consider a rectangular cross-section of a beam. There is the width b and the 
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height h to optimize. Other variables in X come from a model of a general RC 
beam which was presented in [IQ] : the beam is divided to three elements be- 
tween supports with the same diameter of longitudinal reinforcement along the 
top surface and another one along the bottom. The differences are only in num- 
bers of steel reinforcement bars in particular elements. The shear reinforcement 
is designed alike. There are three shear-dimension parts - each of them with 
different spacing of stirrups but the same diameter in the whole element. This 
partitioning reflects the characteristic distribution of internal forces and moments 
in frame structures, where the extremal values usually occur at three points-at 
mid-span and two end joints. The novelty of our approach is the assumption that 
length of parts may attain only the discrete values, in our case corresponding to 
0.025 m precision. The same principle is used for the cross-section dimensions b 
and h. 

2.4 Periodic unit cell construction 

The motivation for this problem comes from the analysis of unidirectional fiber 
composite materials. Such materials consist of a large number of fibers (which 
serve as a reinforcement) embedded in the matrix phase. The goal is to determine 
the overall behavior of such a material system provided that material properties 
of the matrix and fibers are known. It turns out that for this prediction, geo- 
metrical arrangement of fibers must be taken into account. 

Unfortunately, the distribution of fibers in real composite materials is quite 
complex (see Fig. H]). To avoid such an obstacle, we attempt to replace a com- 
plicated microstructure with a large number of fibers by a certain periodic unit 
cell, which resembles the original material sample. More specifically, we describe 
the actual distribution of fibers by a suitable microstructural function and then 
determine the parameters of the periodic unit cell such that the difference be- 
tween the function related to the periodic unit cell and function related to the 
original microstructure is minimized (for detailed discussion see |16j). 

The microstructural function used in the present approach is the second order 
intensity function K{r), which gives the number of further points expected to 
lie within a radial distance r of an arbitrary point divided by the number of 
particles (fibers) per unit area (|14j) 



where /fc(r) is the number of points within a circle with center at the particle k 
and radius r, is the total number of particles (fibers) in the sample and A is 
the sample area. An objective function related to this descriptor can be defined 
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Figure 4: An example of a microstructure of a unidirectional fiber composite. 



where vector = {x^, y^, . . . , , y'^} stands for the position of particle centers 
of the periodic unit cell; x^ and ?/* correspond to x and y coordinates of the i-th 
particle, Hi and H2 are the dimensions of the unit cell (see Fig. (^)), -ft'o(^j) is 
the value of K function corresponding to the original medium calculated in the 
point Tj and is the number of points, in which both functions are evaluated. 
Throughout this study, we assume square periodic unit cell {Hi = H2) and 
determine its dimensions in such a way that the volume fraction of the fiber 
phase in the periodic unit cell is the same as in the original micrograph. An 
example of the objective function is shown in Fig. El^b). 



During last few years, we have developed and tested several evolutionary opti- 
mization methods that are based on both binary /integer and real- valued rep- 
resentation of searched variables. Each of them was primarily applied to one 
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particular optimization problem of the four introduced above. These methods 
are (in order of appearance): 

• Differential evolution, developed by R. Storn and K.Price to solve 
the Chebychev trial polynomial problem. 

• Simplified differential genetic algorithm, developed by authors for research 
on high- dimensional problems 
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• Integer augmented simulated annealing - lASA (a combination of integer 
coded genetic algorithm and simulated annealing); it was primarily applied 
to the reinforced concrete beam layout optimization problem. 

• Real-coded augmented simulated annealing - RASA (a combination of real- 
coded genetic algorithm by Michalewicz [13] and simulated annealing); it 
was developed for solving the periodic unit cell problem. 

3.1 Differential evolution 

The differential evolution was invented as the solution method for the Chebychev 
trial polynomial problem by R. Storn and K. Price in 1996. It operates directly 
on real valued chromosomes and uses the so called differential operator, which 
works with real numbers in natural manner and fulfills the same purpose as the 
cross-over operator in the standard genetic algorithm. 

The differential operator has the sequential character: Let CHi{t) be the i-th 
chromosome of generation t 

CHi{t) = (chait), ch,2{t), chUt)), (7) 

where n is the chromosome length (which means the number of variables of the 
fitness function in the real encoded case). Next, let A be a subset of {1, 2, 
Then for each j G A holds 

chij{t + 1) = chij{t) + Fi {chpj{t) - chqj{t)) 

+ F2{cKestj{t)-chij{t)), (8) 

and for each j ^ A we get 

chij{t + l) = chij{t), (9) 

where chpj and chgj are the j-th coordinates of two randomly chosen chromosomes 
and c/ibestj is the j-th coordinate of the best chromosome in generation t. Fi 
and F2 are then coefficients usually taken from interval (0,1). Fig. [6] shows 
the geometrical meaning of this operator. The method can be understood as 
a stand-alone evolutionary method or it can be taken as a special case of the 
genetic algorithm. The algorithmic scheme is similar to the genetic algorithms 
but it is much simpler: 

1. At the beginning an initial population is randomly created and the fitness 
function value is assigned to each individual. 

2. For each chromosome in the population, its possible replacement is created 
using the differential operator as described above. 

^The determination of A is influenced by the parameter called crossrate {CR), see |15| . 
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Figure 6: The geometric meaning of a certain subtype of the differential operator. 



3. Each chromosome in the population has to be compared with its possible 
replacement and if an improvement occurs, it is replaced. 

4. Steps 2 and 3 are repeated until some stopping criterion is reached. 

As it could be seen, there are certain different features in contrary to the 
standard genetic algorithm, namely: 

• the crossing-over is performed by applying the differential operator ( !8|9|) . 

• the selection operation like the roulette wheel, for example, is not per- 
formed, the individuals that are going to be affected by the differential 
operator, are chosen purely randomly, 

• selection of individuals to survive is simplified to the mentioned fashion: 
each chromosome has its possible replacement and if an improvement oc- 
curs, it is replaced, 

• the mutation operator is not introduced as authors of DE claim that the 
differential operator is able to replace both mutation and uniform crossover 
known from basic GAs. 

Further details together with the source codes of DE can be obtained from 
web page [1]. 
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3.2 Simplified atavistic differential evolution (SADE) 



This method was proposed as an adaptation of the differential evolution in order 
to acquire an algorithm which will be able to solve optimization problems on real 
domains with a high number of variables. This algorithm combines features of 
the differential evolution with classical genetic algorithms. It uses the differential 
operator in the simplified form and an algorithmic scheme similar to the standard 
genetic algorithm. 

The differential operator has been taken from the differential evolution in a 
simplified version for the same purpose as the cross-over is used in the standard 
genetic algorithm. This operator has a sequential fashion: Let (again) CHi{t) 
be the i-ih chromosome in generation t, 

CHi(t) = (chnit), Chi2(t), chinit)), (10) 

where n is the number of variables of the fitness function. Then, the simplified 
differential operator can be written as 

chij(t + 1) = chpj(t) + CR (chqj(t) - chrj(t)) , (11) 

where chpj, chqj and chrj are the j-th coordinates of three randomly chosen 
chromosomes and Ci? is the so called cross-rate. Due to its simplicity this 
operator can be rewritten also in the vector form: 

CHi{t + 1) = CHp{t) + CR{CHg{t) - CHrit)). (12) 

Contrary to the differential evolution, this method uses the algorithmic scheme 
very similar to the standard genetic algorithm: 

1. As the first step, the initial population is generated randomly and the 
fitness function value is assigned to all chromosomes in the population. 

2. Several new chromosomes are created using the mutation operators - the 
mutation and the local mutation (number of them depends on a value a of 
parameter called radioactivity, which gives the mutation probability). 

3. Other new chromosomes are created using the simplified differential op- 
erator as was described above; the whole amount of chromosomes in the 
population doubles. 

4. The fitness function values are assigned to all newly created chromosomes. 

5. The selection operator is applied to the double-sized population, so the 
amount of individuals is decreased to its original value. 

6. Steps 2-5 are repeated until some stopping criterion is reached. 
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Next, we introduce the description of the mentioned operators in detaih 

Mutation: If a certain chromosome CHi[t) is chosen to be mutated, a new 
random chromosome RP is generated and the mutated one CHk{t + 1) is 
computed using the following relation: 

CHk{t + 1) = CHi{t) + MR{RP - CHiit)), (13) 

where MR is a parameter called mutation-rate. 

Local mutation: If a certain chromosome is chosen to be locally mutated, all 
its coordinates have to be altered by a random value from a given (usually 
very small) range. 

Selection: This method uses modified tournament strategy to reduce the pop- 
ulation size: two chromosomes are randomly chosen, compared and the 
worse of them is rejected, so the population size is decreased by 1; this step 
is repeated until the population size reaches its original siz^. 

Detailed description of the SADE method including source codes in C/C++ 
and the tests documentation for the high- dimensional problems can be obtained 
from the article [HI and on the web-page [S]. 

3.3 Real-valued augmented simulated annealing (RASA) 

The augmented simulated annealing method is the combination of two stochas- 
tic optimization techniques - genetic algorithm and simulated annealing. It uses 
basic principles of genetic algorithms ( selection, recombination by genetic op- 
erators ), but controls replacement of parents by the Metropolis criterion (see 
Eq. (1151) ). This increases the robustness of the method, since we allow a worse 
child to replace its parent and thus escape from local minima, which is in contrary 
with DE methods described in Section 3.1. 

The algorithmic scheme of the present implementation is summarized as fol- 
lows. 

1. Randomly generate an initial population and assign a fitness to each indi- 
vidual. Initial temperature is set to Tq = T^ax = T_f racFa^,g and minimal 
temperature is determined as T^m = T_f rac_minFa„g , where Favg is the 
average fitness value of the initial population. 

2. Select an appropriate operator. Each operator is assigned a certain prob- 
ability of selection. 

^Contrary to the classical tournament strategy this approach can ensure that the best 
chromosome will not be lost even if it was not chosen to any tournament. 
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3. Select an appropriate number of individuals (according to the operator) and 
generate possible replacements. To select individuals, we apply normalized 
geometric ranking scheme ([H]): The probability of selection of the i-th 
individual is given by 

p. = q\l - ^r-\ q' = ^ _ (14) 

where q is the probability of selecting the best individual in the population, 
r is the rank of the i-th individual with respect to its fitness, and pop_size 
is the population size. 

4. Apply operators to selected parent (s) to obtain possible replacement (s). 

4a. Look for an individual identical to possible replacement (s) in the popula- 
tion. If such individual(s) exists, no replacement is performed. 

4b. Replace old individual if 

m(0, 1) < e(^(^°'^)-^(^-«/^S (15) 

where F(-) is the fitness of a given individual, Tt is the actual temperature 
and m(-, ■) is a random number with the uniform distribution on a given 
interval. 

5. Steps 2-4 are performed until the number of successfully accepted individu- 
als reaches success_max or selected number of steps reaches counter_max. 

6. Decrease temperature 

Tt+r = TjnultTi. (16) 

If actual temperature T^+i is smaller than T^m, perform reannealing - i.e. 
perform step #1 for one half of the population. 

7. Steps 2-6 are repeated until the termination condition is attained. 



List of operators The following set of real- valued operators, proposed in |T3] , 
was implemented. In the sequel, we will denote L and U as vectors of lower/upper 
bounds on unknown variables, u{a,b) and u[a, 6] as a real or integer random 
variable with the uniform distribution on a closed interval (a, 6). Otherwise we 
use the same notation as employed in Sections 3.1 and 3.2. 



Uniform mutation: Let k = [l,n\ 

c/^.(^+i)=l"^h^^^'^' '[r^- (17) 

^ I chij{t), otherwise, ^ ^ 
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Boundeiry mutation: Let k — u[l, n], p — u{0, 1) and set: 



Lj, if j = k,p < .5 
Uj, ifi = A;,p>.5 (18) 
chij{t), otherwise 



chij(t + 1) = < 

Non-uniform mutation: Let k = [l,n], p = u{Q, 1) and set: 
chij{t + 1) 



' chij{t) + {Lj - chij{t))f, if j = k,p < .5 
chij{t) + {Uj-chij{t))f, if J = k,p>. 5 (19) 
chij (t) , otherwise 



where / = u{0, l)(T't/T'o)^ and b is the shape parameter. 

Multi-non-uniform mutation: Apply non-uniform mutation to all variables 
of CHi. 

Simple crossover: Let k = [l,n] and set: 

ch- (t + l) = I ^^''^^^^ if ^ < A; chji{t), if I <k 

^ ' \ chji{t), otherwise ^ ' [ chii{t), otherwise 

Simple arithmetic crossover: Let k = u[l,n], p = u{0, 1) and set: 

chu{t + l) = I Pchuit) + {I - p)chAtl iU = k 

^ [ chii{t), otherwise ^ ^ 

ch,,{t+l) = lpchAt)+y -p)chu{t), iU = k 

[ chji{t), otherwise ^ ^ 

Whole cirithmetic crossover: Simple arithmetic crossover apphed to all vari- 
ables of CHi and CHj. 

Heuristic crossover: Let p = u{0, 1), j = [l,n] and k = [l,n] such that j ^ k 
and set: 

CHi{t + 1) = CHi{t) +p{CHj{t) - CHk{t)). (22) 

If C Hiit + 1) is not feasible then a new random number p is generated 
until the feasibility condition is met or the maximum number of heuristic 
crossover applications num_heu_max is exceeded. 



3.4 Integer augmented simulated annealing. 

Before presenting the actual optimization procedure we first introduce the map- 
ping between representation and search spaces for individual design variables. 
Consider X = {xi, ^2, . . . , Xn\ as a vector of n variables, integer or real numbers 
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Xi, defined on a closed subset of an appropriate domain Di C Af,7l . Furtlier 
assume tliat eacli variable Xi is represented with some required precision pi, de- 
fined as the smallest unit the number attain. Then, each variable 
be transformed into an integer number ?/j G TV as 

(23) 

where denotes the integer part of XiPi^^. An inverse transformation is 

given by 

Xi = yiPi . (24) 

For instance, the real number 314.159 with precision pi = 0.001 is transformed 
into the integer number 314159. An important aspect of this methodology is that 
the encoded number yi can be treated either as a binary string using bit-based 
operations or as a vector of integer numbers. 

Integer augmented simulated annealing method is based on the same ideas as 
the previously mentioned RASA algorithm. This procedure effectively exploits 
the essentials of GAs (a population of chromosomes, rather then a single point 
in space is optimized) together with the basic concept of simulated annealing 
method guiding the search towards minimal energy states. To avoid well-known 
problems with binary coding the integer coding was used. Together with new 
operators such as differential crossover and a new mutation operator encouraging 
results were obtained. 

The description of the algorithm does not substantially differ from the RASA 
algorithm, but for the sake of completeness all steps are briefly reviewed here. 

1. Initial population consisting of OldSize individuals is created randomly 
and fitnesses are assigned to each individual. Starting and ending temper- 
atures T_niin and T_max are set by the user. 

2. If a real random number p = u{0, 1) is smaller than CrossoverProb the 
crossover is used, otherwise the mutation is applied. This step is repeated 
until the number NewSize of new solutions is obtained. 

3. For each individual in a "new" population one "parent" from "old" part is 
selected. The "old" solution is replaced if 

1 

^) - 1 + e(^(^°ld)-F(/new))/Tt ' (25) 

where F{-), Tt and u(-, ■) have the same meaning as in the previous section. 
Equation ( !25|) ensures the 50% probability of survival when comparing two 
solutions with the same fitness. 

4. Steps 2-3 are performed until the number of successfully accepted individu- 
als reaches SuccessMax or the selected number of steps reaches CounterMax. 



XiPi 
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5. The actual temperature is decreased by 



/ CounterMax 

i = Tt ('^^^\ VTminAtCallsRate * MaxCalls 
\T_ma:xy 



(26) 



where TminAtCallsRate determines a fraction of the maximum allowable 
number of function calls MaxCalls in which the minimum temperature 
T_min will occur. Reannealing step is represented here by setting actual 
temperature Tt+i equal to T_max. 

6. Steps 2-5 are repeated until the termination condition is attained. 

In connection with the notation and principles mentioned in previous sections, 
integer operators within lASA algorithm have the following form: 

Differential crossover: This operator is inspired by the DE. A new individual 
CHjit) is created from three randomly selected solutions CHp{t), CHq{t) 
and CH,.{t) according to 

CHj{t + 1) = CHp{t) + ^(0.0, CR){CHg{t) - CHr{t)). (27) 

Note that all vectors CHi are integer numbers and also that the influence 
of the difference on the right-hand side randomly varies between zero and 
cross-rate CR. 

Mutation: Mutation operator is provided by modifying each variable in CHj{t) 
to 

ch,j{t + 1) = ch,j{t) + N{0, \ ^h^3(t) -^<^K3it) I ^ 1) ^ (28) 

where iV(-, ■) is a random integer number derived from the Gauss normal 
distribution and chpj{t) is the j-th variable of a randomly selected vector 
CHp{t). 



4 Test computations and results 

Each of the methods introduced in the previous section was tested on all pre- 
sented optimization problems. The methodology that has been used for our 
computations is based on the following criteria: 

• For each problem and each method the computation was run 100 times to 
avoid an influence of random circumstances. 

• For all cases, the number of successful runs (which can be traded as a 
probabihty of success or the rehabihty of the method) is presented. 
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• If the number of successful runs is non-zero, the average number of fitness 
calls of all successful runs is also presented. 

Further details of individual function settings and methodology for results 
evaluation can be found in the next subsections. 

4.1 Results for the Chebychev problem 

The computations were performed for the Chebychev problem with a degree of 
the polynomial expression n = 8 (the T8 problem), which corresponds to the 
dimension of the problem 9. The computation was terminated if the algorithm 
reached a value of the objective function smaller then 10~^ or the number of 
function evaluations exceeded 100, 000. Upper bounds on individual coefficients 
were set to 512, while lower bounds were equal to —512. The results of individual 
algorithms are shown in Table [1] and Figure [71 



Method 


lASA 


RASA 


DE 


SADE 


Successful runs 

Average number of fitness calls 


100 
10342 


100 
47151 


100 
25910 


100 
24016 



Table 1: Results for the Chebychev polynomial problem. 



4.2 Results for the type trial function 

Test computations for the type problem were performed for a wide set of prob- 
lem dimensions, ranging from 1 to 200. The upper bound on each variable was 
set to 400, while the lower bound value was —400. For each run, the position 
of the extreme was randomly generated within these bounds and the height of 
the peak uq was generated from the range 0-50. The parameter tq was set to 1. 
The computation was terminated when the value of the objective function was 
found with a precision greater than 10~^. The results are given in the form of 
the growth of computational complexity with respect to the problem dimension. 
For each dimension, the computation was run 100 times and the average number 
of fitness calls was recorded (see Fig. [H] and Table [2]) ■ 

4.3 Results for the reinforced concrete beam layout prob- 
lem 

The basic parameters subjected to optimization were the beam width 6, which 
was assumed to take discrete values between 0.15 m and 0.45 m with the step 
0.025 m and the beam height h ranging from 0.15 m to 0.85 m with the step 
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Number of 




Figure 7: A comparison of results for Chebychev polynomial, reinforced concrete 
beam layout and periodic unit cell problems. 

0.025 m. For each of the three parts of a beam, the diameter and the number of 
longitudinal reinforcing bars located at the bottom and the top of a beam, spac- 
ing and the diameter of stirrups and the length of the corresponding part were 
optimized. Lower bounds were selected for the sake of structural requirements; 
solutions exceeding upper bounds are considered to be irrelevant for the studied 
examples. However, from the optimization point of view, bounds can be easily 
adjusted to any reasonable value. The number of longitudinal bars was restricted 
to the range 0-15, the spacing of stirrups was assumed to vary from 0.05 m to 
0.40 m with the 0.025 m step. The profiles of longitudinal bars were drawn from 
the list of 16 entries while for the stirrups, only 4 diameters were considered. 
This finally results in 18 independent variables. Note that the maximal number 
of longitudinal bars presents only the upper bound on the searched variable; the 
specific restrictions given by Codes of Practice are directly incorporated in the 
objective function. For more details see P,[TO]. The computation was terminated 
if an individual with price smaller than 573.5 CZK was found or the number of 
objective function calls exceeded 1,000,000. Table [3] stores the obtained results 
of different optimization algorithms, see also Figure [71 
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Figure 8: A comparison of results for the type function. 



4.4 Results for the periodic unit cell problem 

Test computations for the periodic unit cell construction were performed for the 
10-fiber unit cell (i.e. the dimension of the problem was 20). The computation 
was terminated if algorithm returned value smaller than 6 x 10~^ or a number of 
function calls exceeded 400, 000. Variables were constrained to the box < Xi < 
Hi 25.8 (see Section (2.4)) . The required numbers of function are stored in 
Table H] and displayed in Figure [71 



5 Conclusions 

Differential Evolution. The Differential Evolution algorithm showed to be 
very efficient and robust for moderate-sized problems, but its performance for 
higher dimensions deteriorated. Moreover, the small number of parameters is 
another advantage of this method. However, the results suggest that the absence 
of mutation-type operator(s) is a weak point of the algorithm. 

Simplified Atavistic Differential Evolution. The SADE algorithm was 
able to solve all problems of our test set with a high reliability and speed. 
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Problem dimension 


T A C A 

lAbA 


"D A C A 




bAUJi/ 


10 


246,120 


13,113 


39,340 


46,956 


OA 


611,760 


74,375 


653,600 


1 '71 r on 

171,539 


ou 


yzD,iuu 


1 S'? 880 
ioO,ooZ 


M / A 


oU4,oZ ( 


100 


2,284,590 


526,492 


N/A 


663,084 


140 


3,192,800 


793,036 


N/A 


948,197 


200 


4,184,200 


1,220,513 


N/A 


1,446,540 



Table 2: Average number of fitness calls for the type-0 function 



Method 


lASA 


RASA 


DE 


SADE 


Successful runs 

Average number of fitness calls 


100 
108732 


100 
131495 


100 
196451 


100 
185819 



Table 3: Results for the reinforced concrete beam layout 



Although it needed larger number of function calls than other methods (see 
Table E]) , the differences are only marginal and do not present any serious disad- 
vantage. Another attractive feature of this method is relatively small number of 
parameters. 

Real- valued Augmented Simulated Annealing. The RASA algorithm was 
successful for all presented problems; the average number of function calls was 
comparable to the other methods. The obvious disadvantage of this algorithm is 
a large number of parameters, which can results in a tedious tuning procedure. 
On the other hand, as follows from the Appendix, only two types of parameter 
settings were necessary - one for the continuous and one for discrete functions. 

Integer Augmented Simulated Annealing. The lASA algorithm was the 
most successful and fastest method on problems with small dimensions. But on 
the problems with larger dimensions and with a higher number of local minima, 
the algorithm suffers from premature convergence and limited precision due to 
integer coding of variables. In addition, initial tuning of individual presents 
another drawback of this method. 

Summary results are given in Table to quantify the overall performance of 
individual methods. Each of the method is ranked primarily with respect to its 
success rate and secondary with respect to the average number of fitness calls. 
The sum then reveals the overall performance of the method. 

Final comments. According to our opinion, several interesting conclusions 
and suggestions can be made from the presented results. Each of them is dis- 
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Method 


lASA 


RASA 


DE 


SADE 


Successful runs 

Average number of fitness calls 


100 
13641 


100 
12919 


100 
93464 


100 
55262 



Table 4: Results for the periodic unit cell problem 



Method 


lASA 


RASA 


DE 


SADE 


Chebychev problem 


1 


4 


3 


2 


Type test function 


3 


1 


4* 


2 


Concrete beam layout 


1 


2 


4 


3 


Periodic unit cell 


3 


1 


4 


2 


S 


8 


8 


14 


9 



Table 5: Overall performance of methods. (* : Not successful for all runs) 



cussed in more detail. 

• The performance and robustness of SADE method was distinguishly better 
than for DE algorithm. This supports an important role of a mutation 
operator(s) in the optimization process. 

• Although algorithms were developed independently, all use some form of 
differential operator. This shows the remarkable performance of this oper- 
ator for both real-valued and discrete optimization problems. 

• The most successful methods, SADE and RASA algorithms, both employ a 
variant of "local mutation". This operator seems to be extremely impor- 
tant for higher- dimensional type-0 functions, where these methods clearly 
outperform the others. 

• Slightly better results of RASA method can be most probably attributed to 
the reannealing/restarting phase of the algorithm (a trivial but efficient tool 
for dealing with local minima) and to the search for an identical individual. 
The procedure for local minima assessment was implemented to SADE 
method (see [71 [S] for results), incorporation into lASA algorithm is under 
development. 

• When comparing methods based on the discrete coding of variables with 
real-encoded ones it becomes clear that for continuous functions the meth- 
ods with the real coding perform better. Nevertheless, after implementing 
new features, like those mentioned before, the performance is expected to 
be similar. On the other hand, the advantage of lASA algorithm is the 
possibility of its use for discrete combinatorial problems like the Traveling 
salesman problem. 



22 



Therefore, from the practical point of view, the SADE method seems to be the 
most flexible alternative due to its simplicity and small number of parameters. 
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Appendix 

See Tables M 



Parameter 


Chebychev, Type 


Beam 


PUC 


pop_size 

Fi = F2 
CR 


10 X dim 

0.85 

1 


11 X dim 

0.85 

0.1 


10 X dim 

0.75 

1 



Table 6: Parameter settings for DE 



Parameter 


Chebychev 


Type 


Beam 


PUC 


pop_size 

CR 

radioactivity 
MR 


lOxdim 

0.44 



0.5 


25xdim 
0.1 
0.05 
0.5 


lOxdim 
0.3 
0.05 
0.5 


lOxdim 
0.2 
0.3 
0.5 



Table 7: Parameter settings for SADE 
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Table 8: Parameter settings for RASA 



Parameter 


Chebychev 


Type 


Beam 


PUC 


OldSize 


80 


900 


180 


200 


NewSize 


5 
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250 


100 


T_max 


10-5 


10-5 


10-^ 


10-1 


T_inin 


10-^ 


10-10 


10-5 


10-5 


SuccessMax 


1000 


1000 


1000 


1000 


Count erMax 


5000 


5000 


5000 


5000 


TminAtCallsRate 


19% 


100% 


25% 


20% 


CrossoverProb 


97% 


92% 


60% 


90% 


CR 


0.5 


0.6 


1.3 


1.0 



Table 9: Parameter settings for IAS A 
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